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RESEARCH ACCOMPLISHMENT IN THE FIRST YEAR 


The research program was planned in four stages over two calendar years 
starting September 28, 1979. The first and second stages would encompass both 
the investigation of free surface properties of the semiconductors in an altra- 
high vacuum (10 -l ° torr) system (UHV) equipped with LEED, Auger electron spectro¬ 
scopy (AES) and ion-sputtering capabilties, and the design and assembly of the 
satellite chamber for device fabrication. These two stages would be performed 
at the same time without interference with each other. The third stage of the 
work would encompass the study of the waveguide surface structure, formulation of 
surface-scattering theory and characterization of scattering loss parameters in 
terms of the surface geometry through LEED and SEM analysis in the same ultrahigh 
vacuum system. The last stage of the work would encompass the fabrication of the 
oxide-free electro-optical devices, the measuring of the optical propagation 
losses and the study of the optical characteristics of such waveguides. 

The first and the second stages of the planned research were completed. 
Currently we are conducting the third stage of the research plan. The satellite 
chamber for device fabrication, as shown in Fig. 1, was designed and assembled 
to the main system. A schematic diagram of the main system is shown in Fig. 2. 

The satellite chamber is indicated as subsystem 2 in Fig. 2. Subsystem 1 in Fig. 2 
is a plasma reactor designed for oxidation research of SaAs and is in no inter¬ 
ference with this research program. A photograph of the completed system is shown 
in Fig. 3. 

The free surface properties of GaAs (110), GaP (110) have been studied by LEED 
and AES. Results of the studies have been presented in two different scientific 
conferences. Abstracts, which summarized the results of the research works, are 
given in Appendices A and B. Abstract of the first paper entitled "Chain Method 
of LEED/MEED Intensity Calculation for Diatomic Surface" is given in Appendix A 


( 







SPECIMEN PLUG 


Fig. 1. Top-view schematic of the satellite system attached to the 
main chamber (left). A: 3orption forepump; B: Vacion UHV pump; 

C: bellows-coupled flange; D: isolation valve; E: six-armed cross 
F: ionization gauge; G: transfer probe magnet actuator. In the 
schematic the transfer probe is shown completely retracted, with the 
specimen plug located at the center of the six-armed cross. 
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f^g. 2. Schematic diagram of the JHV system equipped with LEED, AES, mass 

spectrometer, ion sputtering gun and ellipsometer fpr surface analysis 
Subsystem 1 is the plasma reactor which is currently under constructio 
Subsystem 2 is the MOS device fabrication chamber. constructs 
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LEEO/AES system for Investigation of surface and interface properties of semiconductors 
and ceramic materials. The capital equipment is acquired through the NSE equipment 














































and the second paper entitled "Atomic Structure of GaAs (110) Face" is given in 
Appendix B. Both papers are being prepared for publication. 


The third stage of the research work is currently ongoing. The theoretical 
formulation of surface scattering theory and characterization of the scattering 
loss parameters are currently being developed. The problem of a wave propagating 
along z direction in a rectangular dielectric waveguide with dimensions a and b as 
shown in Fig. 4 was formulated and solved. The wave functions in the waveguide 
were solved using predominantly polarized approximation.* The propagation 
attenuation constants, e.g. the loss parameters, were calculated. Details of the 
theoretical formulation and calculation is given in Appendix C. 

Citations of research results, which are presented in conferences or sub¬ 
mitted for publication, are listed in Appendix D. 


i 


*D. Marcuse, Bell Syst. Tech. J., 48, 3187 (1969). 






Fig. Cross section of a rectangular dielectric waveguide 
with dimensions a and b. ni, n2> n3> n 4 and tv are 
the indixes of refraction of regions 1 to 5, respectively. 
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APPENDIX A 


The following paper was presented in the Conference on Determination of 
Surface Structure by LEED, I3M Thomas J. Watson Research Center, Yorktown Heights, 
N.Y., June 19-20, 1980 

CHAIN METHOD OF LEED/MEED INTENSITY CALCULATION FOR DIATOMIC SURFACES 

In low energy electron diffraction (LEED), crystal surface is probed by bom¬ 
barding electrons of energy E <_ 200 eV at normal incidence to the surface and analy¬ 
zing the intensity of those elastically scattered electrons which are back reflected 
from the surface, whereas in medium energy electron diffraction (MEED) due to an in¬ 
crease in the energy of incident electrons (200 eV < E <_ 5 keV) the surface sensiti¬ 
vity of the technique is maintained by taking an oblique angle of incidence. Inter¬ 
pretation of LEED/MEED intensity spectra requires detailed calculations of the 
diffracted intensities for a series of trial models. For such a technique to work, 
a rapid and accurate method of calculation is needed. 

Layer-KKR method is conventionally employed in the calculation of LEED intensity 
spectra. In this formalism, the crystal is divided into a number of layers parallel 
to the crystal surface and scattering calculations are split into two parts, intra¬ 
layer and interlayer. In the intralayer scattering calculation, due to the assumed 
spherical symmetry of the atomic potential, an angular momentum representation is 
used. This involves large matrices and their inverses and becomes cumbersome with 
increasing energy and more complex surfaces. In such a situation we propose to use 
the chain method of intralayer multiple scattering calculation. 

In the chain method the two dimensional intralayer multiple scattering calcu¬ 
lations are further divided into two one dimensional steps, scattering within a 
chain of atoms and between the chains. At each stage of calculation the electron 
wave function is represented in terms of an appropriate set of basis functions. The 
scattering by an atom, a chain of atoms and a layer of chains is expressed in terms 
of spherical, cylindrical and plane wave representation respectively. As we move from 
one stage to another, transformation from one basis set to another is carried out. 

Due to its one dimensional lattice summations the chain scattering formalism has 
many computational advantages over the layer-KKR method both for normal incidence 
LEED and off normal incidence MEED. 




APPENDIX B 


The following paper is submitted for presentation in the Eighth Annual Con¬ 
ference on the Physics of Compound Semiconductor Interfaces, Williamsburg, Virginia, 
January 27-29, 1981. 


ATOMIC STRUCTURE OF GaP (110) FACE 
ABSTRACT 

Low energy electron diffraction (LEED) intensities have been measured for the 
(110) face of GaP and analyzed using a dynamical multiple-scattering model of the 
diffraction process. The intralayer multiple scattering is treated exactly, while 
for the interlayer multiple scattering, the renormalized-forward-scattering method 
is used. Comparison of the calculated and observed LEED intensities suggests that 
both the Ga and the P atoms on the (110) face may exhibit a contracted outermost 
layer spacing. The surface layer is compressed by about 5% such that the top layer 
spacing is reduced by 0.1 ± 0.02A. The rippled geometry of surface reconstruction 
is not clearly observed. This indicates that the GaP (110) surface atomic struc¬ 
ture is different to that of GaAs (110). By comparing the structures and properties 
of the GaP (110) face and the GaAs (110) face, it is concluded that the GaP (110) 
face is relatively unstable and reactive. 









APPENDIX C 


The problem of a wave propagating along z-direction is a rectangular dielectric 
waveguide with dimension a and b as shown in Figure 4 is formulated and solved. An 
exact analytical treatment of rectangular waveguide is practically impossible. 
Therefore the approximated analytical approach developed by Marcatili (a) is followed 
in solving the problem. 

Assuming a traveling wave propagating in z-direction as shown in Fig. 4, 

Maxwell equations 


7XH = V 2 f 


’ u . I 


can be solved to give 
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where n is the index of refraction of the dielectric; k is the wavevector in 
free space; u> is the angular frequency of the wave and 3 is the propagation 
constant in z-direction. 

Substituting (2) and (4) into (5), and (1) and (3) into (6), we have 


3 2 E 3 2 E 
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Two different types of modes can be supported in the waveguide v , e.g., 


: polarized predominately in x-direction and 



E* : polarized predominately in y-directi on, 

r i 

where p and q are positive integers indicating the modes of the wave propagating 
in the waveguide. 

E* MODE 

pq 

In region 1 with refraction index n^ as shown in Fig. 4, the waveequation (7) 
can be solved to give 


E z (x,y) = A Cos k x (x+c) Cos k y (y+n) 


with = 0, where £ and n are the phase factors of E z> 

Other components of the wavefunction can be obtained by substituting (9a) 
into eqs. (1), (3) and (4): 


H - - fiX- (% VZ (£*) Sin k (x + C) Sin k w (c + n) 


3 v u ' v k ' 
H o x 


(n?k 2 -k 2 ) 

E x = iA Tk-” Sin k x + ^ Cos k y + ^ 


E y = —g- 2 - Cos k x (x + c) Sin k y (y + n) 


Substituting (9a) into (7) we have 


n|k 2 -0 2 = k 2 + k 2 
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For small incident angle, 3»(k 2 + k 2 ), eq. (10) can be approximated to be 

x y 

3 = n.jk Cos9 = n^k for 9<5°. 

The wavefunctions in regions 2, 3, 4 and 5 with indices of refraction n^, 
n^, and n^, respectively, as shown in Fig. 4, are solved. When boundary con¬ 
ditions are matched between region 1 and the neighboring regions, i.e., regions 2, 
3, 4 and 5, we obtain the four transdental equations: 


tan k x a = n 2 k x (n 2 y s + n 2 y 3 )/(n 2 n 2 kj - n*Y 3 Y 5 ) 

tan kC = (^ (k x /y 5 ). 

tan k y b = k y (y 2 + Y 4 )/(k y 2 - Y 2 Y 4 ) . 

tan k y n = -Y 4 /k y . 

where y 2 = [(n 2 - n^k 2 - k 2 ] 1/2 

v 3 = C(nf - n|)k* - 


(Ha) 

(lib) 

02a) 

(12b) 


Y 4 a C( n i - ] V2 

Y 5 - C<"? - - k J ] 1/2 

A computer program is currently being developed using Newton's method to solve 

equations (11) and (12) to obtain k , k , 6, and the power attenuation constant a. 

x y 

REFERENCES: 

(a) Marcatili, E. A. J., Bell Syst., Tech. J. 48, 2071-2102 (1969). 

(b) D. Marcuse, "Theory of Dielectric Optical Waveguide", Academic Press, New 
York, 1974. 

(c) I. P. Kaminow, W. L. Mammel, and H. P. Weber, Applied Optics, Vol 13, No. 2/ 
Feb. 1974. 
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APPENDIX D 


Publication citations and presentations supported or partially supported by 
this research grant: 

1. "Decomposition of Aluminum oxide by Electron Bombardment." B. W. Lee and 
J. M. Kuo, 8. Am. Phys. S ., 25 (3), 238, 1980. 

2. "Study of MIS Polycrystalline Silicon Solar Cell Using Auger Electron 
Spectroscopy," J. M. Kuo, 8. W. Lee, B. Lalevic and W. A. Anderson, B. Am. Phys. S ., 
25 (3), 409, 1980. 

3. "Chain Method of LEED/MEED Intensity Calculation for Diatomic Surfaces," 

N. Masud, C. G. Kinniburgh, D. J. Titterington. Presented in Conference on Deter¬ 
mination of Surface Structure by LEED, IBM, T. J. Watson Research Center, 

Yorktown Hights, NY, June, 1980. 

4. "Study of MIS Silicon Cell by ESCA and AES," Y. S. Wang, H. J. Yu, C. C. Hsu 
B. W. Lee and W. A. Anderson, presented in the 27th National Symposium of American 
Vacuum Society, Detroit, October, 1980. Submitted to J. Vac. Sci. Technol for 
publication. 

5. "Stability of MIS Silicon Solar Cell," B. W. Lee, J. M. Kuo, B. Lalevic and 
W. A. Anderson, submitted to J. Vac. Sci. Technol . for publication. 

6. "Atomic Structure of GaP(110) Face," B. W. Lee, R. K. Ni and N. Masud, 
submitted for presentation in the Eighth Annual Conference on the Physics of 
Compound Semiconductor Interfaces, Williamsburg, Virginia, Jan. 1981. 
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A NONLINEAR MAXIMUM ENTROPY METHOD FOR SPECTRAL ESTIMATION 

Summary 

An intensive research study is made of the nonlinear maximum 
entropy spectral analysis method proposed by P. F. Fougere of 
AFGL. The research not only provides a much better understanding 
of the properties of the method including spectral resolution, 
convergence, etc. but also shows that the method can be performed 
effectively with the PDP 11/45 minicomputer. The successful 
implementation of the method is described in the program listings 
in Appendix I. Extensive computer results are presented or. 
various data. These results clearly confirm that the nonlinear 
method is superior to the Burg's maximum entropy spectral analysis. 

Preliminary results on the multichannel (multivariate) maximum 
entropy spectral analysis, the two-dimensional maximum entropy 
spectral analysis, and computer graphics for the spectral display 
are also presented. An extensive bibliography of the maximum 
entropy spectral analysis is given in Appendix III. 
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I. Review of Research Progress 

This research is concerned with the nonlinear maximum entropy 
spectral analysis (MESA) method proposed by Dr. Paul Fougere of 
AFGL. As verified in this research work performed at the 
PDF 11/45 minicomputer, the method not only provides a much better 
spectral resolution than the Burg's method but also removes the 
line-splitting and frequency shifting phenomena for sinusoidal 
signals, as experienced in the Burg's method. By using the double 
precision, the minicomputer results are reasonably close to those 
obtained at the CDC 6600 computer at AFGL. 

The final computer program developed for the nonlinear 
complex signal maximum entropy spectral analysis is shown in 
Appendix IA. The program follows the mathematical development of 
Fougere [*l] but is quite different from the original computer pro¬ 
gram developed by Pr. Fougere. Appendix IB is the computer 
program for Burg's complex signal maximum entropy spectral analysis. 
By using the computer programs and the two-channel radar aata as 
shown in Fig. 1 (see also["2] ), the spectrum of the nonlinear 
method is shown in Fig. 2a (linear plot) and Fig. 2b (logarithmic 
plot) for 10 filter weights. The 3urg's result is shown in Fig. 3a 
(linear plot) and Fig. 3^> (logarithmic plot), also for 10 filter 
weights. More detailed tabulation of the major frequency compo¬ 
nents outside the clutter bandwidth (-0.l6?fs, +0.l67fs) is given 
in Fig. 4. The nonlinear method which is very close to the true 
answer clearly is much better than the Burg's method. 30 iterations 
are used in the nonlinear method which appears to be an optimum 
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number. The optimum filter weight is around 10 or 11 as 
determined by Fig. 5 which shows the finear prediction error (FEZ) 
according to the Akaike criterion. The solid curve is the lower 
bound and the dashed curve is the upper bound as tabulated in 
Fig. 6. The upper bound was originally proposed in the research 
proposal [ 2 ] and it seems to be better than the lower bound due to 
Akaike. 

Major documentations already made which describe the research 
progress are as follows*. 

1. C.-H. Chen, J. Chen, and C. Yen, "A minicomputer implementation 
of Fougere's maximum entropy spectral analysis method," Techni¬ 
cal report prepared for the Mini-grant, August 20 1980. This 
report has detailed results of comparison between Burg's and 
the nonlinear methods for sinewave and sunspot data. 

2. C. H. Chen, "Spectral resolution of Fougere's maximum entropy 
spectral analysis," to be published in the Proceedings of IEEE, 
June 1981. This journal paper based on the work performed 
under the Mini-grant provides a good comparison between the 
nonlinear and Burg's methods, and the other method for complex 
sinusoids. The Cramer-Rao bound is used as a reference. 

3 . C. H. Chen and C. Yen, "Note on computer graphics for maximum 
entropy spectral analysis," Technical report prepared for the 
Mini-grant, March 23> 1981. This report provides a three- 
dimensional spectral display of sinewaves for both nonlinear 
and Burg's methods. 
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A number of important results are included in the new research 
proposal submitted to AFOSR in December 1980. Appendix III 
provides an extensive list of references on the maximum entropy 
spectral analysis. The following sections describe some new 
research areas with preliminary results. 

II. Multichannel (multivariate) Maximum Entropy Spectral Analysis 

A mathematical presentation of this topic is given in Appendix 

II. Several computer programs for multichannel maximum entropy 
spectral analysis were provided by Dr. Fougere. The following 
results are based on the time series of sunspot numbers, northern 
light activity, and earthquake activity by using the third multi¬ 
channel program. The data are tabulated in [3]- 

Fig. 7a is the first channel (sunspot number) auto-spectrum 
with linear (left) and logarithmic (right) scales, and 16 lags. 

Fig. 7b is the second channel (northern light activity) auto- 
spectrum with linear(left) and logarithmic (right) scales, and 
16 lags. 

Fig. 7c is the third channel (eartnquake activity) auto¬ 
spectrum v/ith linear (left) and logarithmic (right) scales and 
16 lags. 

Fig. 7d is the cross-spectrum between channels 1 and 1 with 
real part (left) and imaginary part (right), and 16 lags. 

Fig. 7e is the cross-spectrum between channels 2 arid 3 v/ith 
real part (left) and imaginary part (righ ), and 16 lags. 

Fig. 7: is the cross-spectrum between channels 1 and 3 v/ith 
real part (left) and imaginary part (right), and 16 lags. 

Dy way of verification, it is interesting to note that the 
spectral peax for the sunspot numbers is determined accurately. 

III. Two-Dimensional Maximum Entropy Spectral Analysis 

be consider a very simple separable case in the two-dimensional 
spectral analysis. The signal considered is sin(2Fx) sin(2ITy). 

In t.uis case the power spectrum is the product of the power spectra 
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of sin(27Tx) and sin(2Hy). In each spatial dimension, the power 
spectrum can be determined by using the Fourier, nonlinear and 
Burg's methods. The two-dimensional spectra are shown in 
Figs. 8, 9 and 10 respectively, based on the Fourier, nonlinear, 
and Burg's methods. The nonlinear method clearly is much better. 
Extension of the above procedure to a more general two-dimensional 
spectral analysis is not possible. Although some two-dimensional 
maximum entropy spectral analysis work has been reported (see 
Appendix III), the success is very limited. Further research is 
much needed. 

IV. Conclusions and Recommendations 

The nonlinear maximum entropy spectral estimation method 
proposed by P. F. Fougere has provided superior spectral estima¬ 
tion over the Burg's method in a number of data considered. The 
computational requirement of the nonlinear method is, however, 
significantly higher. Typical number of iterations needed is 20 
to 30. The use of minicomputer has not created significant 
computational problem as predicted. It is our strong belief that 
the nonlinear method will become very popular in high resolution 
spectral analysis for a wide range of applications in geophysics, 
sonar, radar areas, etc. 

Recommended further studies include the multichannel 
(multivariate) maximum entropy spectral analysis, the two- 
dimensional nonlinear maximum entropy spectral analysis, computer 
graphics for the spectral display, signal decomposition, and 






signal prediction and extrapolation. 
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No. of weight 

Log. of 
N+(M+1) „ 

Log. of 

N+2(M+1) 

N-(M+1) F m 

N-2CM+1) P m 

1 

-0.279 

-0.224 

2 

-0.463 

-0.381 

3 

-0.445 

-0.332 

4 

-0.562 

-0.418 

5 

-0.674 

-0.496 

6 

-1.080 

-0.866 

7 

-1.539 

-1.283 

8 

-2.063 

-1.762 

9 

-2.587 

-2.006 

10 

-2.597 

-2.175 

11 

-2.678 

-2.176 

12 

-2.686 

-2.075 

13 

-2.738 

-1.971 

14 

-2.815 

-1.764 

15 

-2.839 

- 

16 

-2.824 

- 

17 

-2.833 

- 

18 

- 2.810 

- 

19 

-2.799 

' 


Lower Bound Upper Bound 


Note: N is the number of data points in each channel (32 in this case) 
M is the number of filter weights 
P m is the error power 


Figure 6 
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Appondix IE Burg's Complex Signal. MESA Ccmputer Prcgrarn 


OKI KAN 

006. 

i 3 1.10 o 3 4 1'3-14 A K - 8 i r A13 F 


c 

IMA) iF : CMPBUO C 11. 1 U OATF : 14-AUG-80 

>001 


COMMON IN)Hr X 

■ >002 


REAL F(32),PbF(32),PFR(32),FNkW(32>, CNFN(32) 

1003 


Hh:.Al Vi (70 ) i Vi Vi Vi ( 70 ) > OViVill ( 70 ) < HAi. F ( 1 f . 1 ) ■ H ( 7) 

•004 


Kt-.Al. SI- (.313) , 8(1)13). YX( 1O00) < CAL< ' 00<>) 

i >005 


UU! iPLKX h/ Phhi F L K j Vi Vi 1.7 Vi Vi ij i'll 0; FU-il.F* H 

•003 


1<!K I 1 F ( 6< 3 ) 

:>007 


FORMAT ( 1 X, ' I/p NO 11. I.G. N EXPAND, )j'i #-*•*•» 

>003 


KKAO(6, lO/NUII, 1 0, Si !. EXPAND, DT 

..<009 

10 

FORMAT (213, 3F :■) 

**>0 ). 0 


Wh]Tr.<6, 15) 

.'01 A 

13 

FOR 110 1 ( 1 >., '*•#**•# !./-• IS 1 AN 1 , IS ! UK, INC, | in j l , 12 

0012 


RK0.0 (6, 20) 1.81 ART , IS I'UP, INC, NO' M , I 2, I.K 

i013 

20 

K OK'i'iA'f ( 6 _• 3 > 

.'01'I 


WRIT E(3, 33) NUI1, LG, I STAR 1 , Is COR, NIC, NOLL , 12, SN 

•015 

\iL 'I' 

FORMAL( 1H ' , 3X, /17, '21? 3. 

0013 


CALL EFLI. 

O01 7 


N : - 4011*7 

•..1013 


DLK INL K 1 LE 2 ( 3, NI, IJ J. N) H- X ) 

0019 


CALL INPL T ( K , HUM S N, K IMF k, GNE.N ) 

.>020 


WRI I F. (5, 30) NUN 

0021 

30 

H* VJFI'jA T (c >X i F DA 1 A 1 , U. v.*F FF - > - , [-"■) 

0022 


CALI. WR ( NUN, F , 5) 

0023 


IKr-8. *ATflN( 1. > 

•0024 


CD! L OKU 

0025 

50 

CONr1NUK 

0026 


DFN-FLOAi (NCU l-1 ) ttEXPAIID 

(.) 02 7 


LFLU..' - IN 1 ( OFN/KXF'AND -1. > 

0028 


l.i- • 60 KF- 1, LFUI 

C 029 


SF ( K F ) —F L UA 1 ( Li.1. ) /) (E.N—. 3 

0030 

60 

SF (KF ) •-! PI ■»SF ( KF') 

0031 


DO 63 I-I.IMLHI 

’ 0032 


FI'EM ( I >-KFAL(F ( I ) ) 

• 0033 


GNK1J ( I '■ “ATI4A0 (E ( I ) > 

0034 

4 5 

C( IN r ] N JF. 

! 0035 


DU 70 I-5MLH1 

0036. 

70 

YX < I ) : =F LG AT ( I : 

i 003/ 


LAI L NEWE'AG 

0038 


L 01 L LERI Ol (YX, KNEW, NUN, 50, 1010, 60, 700' 

0039 


CAI L L IG1 ( YX, CAL, 348, 3) 

(1040 


CAl 1 Eh 1 1 

1 004 1 


OKI T E ( 6, Cl ) 

0043 

y v:£ 

K OKI' AT ( 20X , -UHHMHt I / p CHLI'INFL i WAVE *•**#* 

0043 


REA) 1 < 6, 999 ) NX 

004 4 

yyy 

KUSMA1(13; 

004 5 


i KIL L NEW RAN 

0046 


CALL KEPI NT (YX, C-IMR- • NON !;(>, 1010, 3C. , /CO) 

004 / 


CALL FIG1 (YX, LAL. -48.3) 

0048 


CALL EFI 1 

0049 


NP L I E ( 6, 9 -5 ) 

0050 

9 ? ,’ 

FOkMA r (20X, ''*■»**•*■ 1 /-■ CHANNkl I t WAVE *•**■»■* •' 

0051 


READ (6, 99 , -')HX 

0052 


F - 0. 

005 3 


DO SO I - 1 , MUE 

i 0054 

8C> 

P-.-M-F ( I ) *>: OiM.JG < F ( I ) ) 


/ 
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ORTRAN 

VOo. IS 

00 SO. 14 2b—MAK-8' 

■ :»55 


P~E’/ELC'Al (NUh) 

•056 


SAVE=P 

057 


WRITE(5,£b>P 

‘058 

X ”i 

E UK El AT (21 X, •" P - ■ , E 1 b 7 ) 

i.i 5,9 


LA 1 -NO El -1 

••:)/-,() 


DU 150 NIL-i, LA J 

CM 


NNN-NN+1 

.'062 


CALL.. 8H‘H_L‘ ( MUM/ N : \l > H*H , H'HKi H» G# GGGi GGGi v i> H) 



P-( 1. -G(I'll'.! >■wCON-JG ( C- ( NK > ) > *P 

'064 


HPl.P-P 

065 

;i 00 

IF ( NNIM u 1 .1 S i (ip) i iij TO 1 bO 

v,6/-. 


IE < MMN. L I 1ST AP I ) GO TO ISO 

'067 


IPUIOLK (NUN--1ST ART '•, INC). ME. 0)00 ' 0 V 50 

06s 


L AI.L ’r*PhL*(NN» i.L Li ’"'Hi KHUI i r > D f ) 

o/.v 


CAM BEI...I 

0/0 


Col 1... BE LI 

071 


READ(6,99 9>NX 

>072 


CAI 1. NEWPAG 

•0/3 


CALL KBPI.GRYX, S, NulJT, 50, 1010, 50, '00) 

>074 


WRITE(6, 170>MNM 

‘07 b 


WRITE(5, 1 "(0 ) Nlv 

•076 

3 70 

H UK HA 1 ( 80 X * ‘ ->r & -i(- |vQ. OH HJHIGHI - '* • 13. 

• 077 


CAI L UR(NN, G, 5> 

•078 


CAI L PIG) (YX, CAL, 12, 7) 

>079 


IE (IP. KQ. :l. ) CAI. 1.. PR IT ( NGU i • DT , S, Y > • 

>080 


IE (LG. LE 0 ) GO TO 1 40 

>081 


R-- AD (/., 9'- , 9 ) NX 

>082 


CAI.! AUKS, NOLO ) 

.‘083 


CALL MENRAG 

>084 


CALL WR(NN,G, 5) 

>085 


CAM. KBPI..01 (YX, S, NulJT, l;0, .1010, 50, 700 • 

>086 


WRITE (6, 170)NIMM 

>087 


CALI FIG1 (YX, CAL, 12, 1 > 

.•088 

140 

E-HELP 

>089 

150 

CON 1 INI IE 

>090 


CAI L BE.LI 

>091 


CALL EXIT 

>092 


ENU 


KOI.IT INK 

S CALI Eij: 


BELL , 

INPUT , UR , AIAM , EI.OAl , Jl-.ll 


AI MAG , 

MENRAG, KBRLOT, E 101 , CON..1C- , BREL 


SPEC , 

PR IT , Al G , EXIT 


CRT] iINS 

= /0P: 2 


bi .o cl: 

LEMO 1 H 


MAIN. 

81S5 (037762)* 


. 

2 (000004) 


♦♦COMRILER-CORE#* 

PHASE USED EPEE 

DECLARATIVES OOS-SA 1451? 
E X El '-LIT ABI.hS 0126 3 1411 b 
ASSEMBLY 01/74 18044 


W-V-.v 


/ 


I . 


1 

i 

i 


i 

i 

i 

■i 














i 
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I 

> ' 'RTKAN VOo. 


■ >01 
\ >02 
i i03 
>04 
■ >0b 
>06 
>07 
• •: >08 
■o 00 
010 
•oil 
>01 2 
•.<13 

•Ml 4 10 

•015 

016 20 

•01V 

• OJ 

•019 30 

■020 40 

021 
<022 


00 >1 4 4 25-l / -AR-31 

SUBROUT INK INKUT(X,Ml,SN, Y, L) 

COMMON INDEX 
COMKLb X XCi) 

REAL Y(N ■ ), X (N ) 

iuio 

READ (2 INDKX )Y 
I NDb X= 2 

Rh.AD(2 " INDEX ) Z 
DO 0=1,01 

z —z •: o > 

CO NT INI 111 
DO 10 I 1 I' 1 1 

X( I )= CMKLX (Y (I >, Z (I ) ) +SN*CMM. X (GAO'::S( 
CON 1INUb 
NR 1 ! h: (' j i 20 ) 

F OR HAT ( / "AX, CHANNEL I ", ' <X, •' CHANIMF i > 

D>j 40 1 -1 , Mi 

UK 1 I E (’: • SO > Y ( I >, Z i I ) 

K OKI'I A1 ('."'EX K V. Ei b X, FV. A) 

CCNI INOE 
RHIURN 
KND 


ROM I I NFS CALLkB: 

CMKLX , GAUSS 

OKI IONS -/ OK: 2 

61 OCT . Lb NO l H 

INKUT 2 XT (001042) ■>(• 

. *$$•!>. 2 (00000 £ ) 


^COMPILER — 
F'KASb 

DECT AKATlVbS 
b XF.CU I AGLb.S 
ASS EM 61 .Y 


- C0RbOO> 

USED b Rbb 
0062: v 14 7'S ', 
00736 143V2 
01218 13800 


1 

l 



/ 


KA> <b 


1. i (..• < i 6AL 

I ' < 
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'RTKAN 

VO/-. IS 

00: 33: 11 2b-NAP 

oOl 


SUB ROUT I. nr NR; MPP 0, IC) 

..>02 


Rh AL G ( 1 ) 

003 


COMF1 p:x g 

- 04 


DO 10 I ~ 1 , NN 

00 b 


A-RE At . (G< I ) ) 

>06 


R ; "A I MAG (0 ( I ) > 

00 7 


[PCI C. PIO. 5 ) WRI 1 P" (5, 20 ) A, B 

008 


IP ( 1C. POT. 6)UKIYE(o. 20)A. B 

009 

20 

FORMA f (28X. ( ", P 15. 7, ", P lb. 7, " ) - 

01 0 

10 

CONTINUE 

oil 


RETURN 

012 


e m 


ROUTINE 

S CAL.LEU: 


REAL , 

AI MAG 


AP I IONS 

-/OF: 2 


Bl OAK 

LENG'I H 


WR 

] bO (000454 )-> 


»‘#COMKILEK - CORM-» 

PH Ash USED h Rh £ 

11KCLARAT ) V- S 0 <>7 0> 1 46 7 6 
h Xh CU ! API HS 00783 1 4b9b 
ASSEMBLY 0109/ 13921 


FORTRAN 

V06. 5 3 

00:32:31 

ooo :i 


P UUCTI ON GAUSS ( XBAK, SI Gl 1A > 

(>002 


ROC r * 1 / A Uhl ( S -lATAN ( 1. )) 

OOO/- 

10 

corn i hup 

i)004 


X 0 •)■• ( RAN ((>i 1 pi) 

0005 


Y -RUOTS-tXP ( —X#X/2 ) 

0006 


YTRY'-ROCi'l *RAN< 0. 1 ) 

0007 


IPCYTRY tTf Y) GO 10 10 

000:3 


G A: i s s=X » SIG M A ■+ X B A R 

000V 


RP TURN 

0010 


PND 


Riil.If INKS LAI..LED: 

SORT , ATAN .. RAN , tXP 


MR riLNS = 'f’P 2 


Bl III • I PNG I H 

■ »AUSS 1 'X (00'>R 2S ) •«• 

**1 i.MPILEK-CLRP-i 1 * 

PHASE l ISPD PREP 

DP i. I .AKA i i VI- S * >06 '9 \ 4 T <k 
EXRi lllABLPS 00" S 3 1459 b 
ASSEMBLY 0103 * IS-P'A'j 












ORTKAN 

006 13 

00:32:56 25-NAK-81 PAG*- 1 

001 


SUBROUTINE BP EC ( NUI-I. NN, PPP, PPR, P, 0, GAG, nGGfi, H • 

i002 


DIMENSI ON H ( 1 ) i PP P ('!.)< P'P P: ( 1 ) , P ( ) » ■ i ( 1 > < U> jG ( '• ) . • •>> •• .41 ( 1 

'003 


COt* PI h.X U# PKh i FFK. H. OViCl v i» F# F'. U FM. * 

•004 


XP.RC»=<0. , 0. > 

• 005 


N--NM— 1 

•006 


IP (N. NPI. 0)00 TO 20 

>007 


DO 10 J~1,NUM 

>008 


PPP (J)>»XP.KO 

•009 

10 

PBK(J>=XhRO 

>010 

80 

SN-XP.RO 

•Oil 


SD -XPRO 

•012 


J J- NON—N— 1 

• 018 


DO 30 J~1,JJ 
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Appendix II 

On Multichannel (Multivariate) Maximum Entropy Spectral Analysis 


1. Introduction 

The univariate maximum entropy spectral analysis has now been well developed 
and applied to many defense research areas such as radar, sonar and geophysics. There 
has been some work done to extend the maximum entropy spectral analysis to multivariate 
case. Whittle [1] and Robinson [2][3] generalized the Levinson-Durbin recursion to 
the multivariate case by fitting both foward and backward autoregressions in a 
stepwise fashion. In this thesis. Burg [A] has mentioned about the multichannel case. 
However the computer programs for both multichannel and multivariate maximum entropy 
spectral analysis were only recently developed successfully. Morf et.al [5] developed 
an algorithm for direct estimation of the normalized reflection coefficients from the 
observed data for maximum entropy spectral analysis. They also compared the spectral 
estimation with the methods of Jones [6], Nuttall [7] and Strand [8], which are more 
of a direct extension of Burg's work to the multichannel (multivariate) case. Burg's 
algorithm does not generalize directly since the forward and backward autoregression 
matrices are not the same in the multivariate case, and the forward and backward one- 
step prediction error covariance matrices are different, although they have the same 
determinant. In this report, the programs developed by Strand and Jones are applied 
to real multichannel data and imagery data in addition to a set of test signals. 

The merits of these methods are closely examined. In spite of programming complexity 
the multichannel and multivariate maximum entropy spectral analysis will have increased 
application as the real data are almost always gathered in several channels. Data 
from several channels form a vector for multivariate study. 

2. Brief Mathematical Analysis 

Let x r x 2 ,...,x n denote n zero mean vectors of dimension d each. The sample 
estimate of covariance sequence for lag j is 


/ 
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i n ; j 

R. = — l x x' 
J n t ^ 1 t+j t 


( 1 ) 


where the prime denotes the transposed vector. The forward and backward predicting 
autoregressions of order p are given, respectively, as 

\-k 
(P) 


x (f) = l A J^x 
C k=l k 


-Cb) = 


^ = l B 
C k-1 k 


( 2 ) 


t+k 


where A R P ^ and B R P ^ are d x d matrices, and can be determined recursively [6] by 


making use of the estimated covariance matrix in Eq. (1). The recursion starts with 

s ( f ) = g (t) = R 

0 0 0 

The one-step forward and backward prediction error covariance matrices are 


(3) 


S (f) = (i _ A ( P ) B ( P ) )S (f) 1 
p P P p-1 

S (b) - (I - BWA ( >V b > 

P P P * P~1 

The forward and backward residuals are, respectively, 

„(p) 


(A) 


x - l A< p) x 
t , k t-k 
k=l 


t = p+1,...,n 


(5) 


3 t P) *t J- B k P)x t+k’ C p 

k=l 

The recursive equations are then given by 


-? 1 

= e £ p-1) 

_ A (P) 

P 

gCp-D 

e t- P • 

t = p+l.n 

(6) 

6 (p) 

- 

- B (p) 

.Jr 1 ’. 

t = 1,...,n-p 


t 


P 

t+p 



least 

squares 

estimates for the 

forward and backward 

autoregression matrices 

A (p) 

= uv- 1 




(7) 

P 




B (P) 

= u’vT 1 






where U is the sum of cross products of forward and backward residuals at lag p, 

u-7 W 


t=l 


’t+p 


■m -‘♦t - - v 


- 
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and V and W are estimates of (n - p)S^ > ( n _ p)s ( '^ respectively, 


n-p 


V = I 6^ 

t = l 


(p-D /g(p -1 )\ i 
t 1 


(9) 


w- Y e'r 1 ’ <Jr l v 

t+p t+p 


( 10 ) 


Although the forward and backward autoregression matrices and the prediction error 
covariance matrices are different, the multivariate spectra should be identical 
when calculated from the forward and backward fits by 


S(f) = h[A(f)] _1 S< f) [A*(f)] 

or by 

S(f) = h[B(f)] _1 S^ b) [B*(f)] 

where 

A(f) = I - l A< p) e 2lTikhf 
k-1 k 


B(f) - I 



-27rikhf 

e 


h is the sampling period and * denotes complex conjugate transpose. 

The above approach based on the work of Jones [6] does not guarantee stability 
and does not generally produce a non-negative definite spectrum as has been pointed 
out by Nuttall [7], Subsequently, Nuttall [7] and Strand [8] applied a weighted 
arithmetic mean error criterion in order to provide model stability and to ensure 
positive definite stationary spectra. Another procedure suggested by Morf, et.al. 
[5] that also meets the spectral requirements is to compute the spectrum from the 
normalized reflection coefficient matrix p. To obtain this matrix, W and V are 
factored by using Cholesky decomposition into the product of lower triangular 
matrices times their transposes. A new recursive procedure for and by 

using p in place of Eqs. (A) can then be obtained. Other recursive algorithm has 
been proposed [9] for the solution of the normal equations for both single and 



/ 



multichannel data. 
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Appendix III 


A BIBLIOGRAPHY ON MAXIMUM ENTROPY SPECTRAL 
ANALYSIS AND RELATED TECHNIQUES 


I. Introduction 

Recently there has been strong research interest on high 
resolution spectral analysis techniques. This is an important 
area of defense research because of the numerous applications to 
radar, sonar, and geophysical areas of defense interest. An 
excellent publication is the Proceedings of the 1978 and 1979 RA^C 
Spectrum Estimation Workshop. Maximum entropy spectral analysis 
is one of a number of high resolution spectral analysis techniques. 
The impact of the Burg's maximum entropy spectral analysis method 
is far more significant than the technique itself. Thus in this 
report we present not only the bibliography of the maximum entropy 
methods in one and two spatial dimensions but also a number of 
related methods of high resolution spectral analysis. Cne common 
assumption with all these methods is that the data record is short 
and thus the conventional fast Fouriter transform method of spec¬ 
tral analysis is not suitable. Frobably because of the snort length 
record, the maximum entropy spectral computation is fairly sensitive 
c the presence of noise. In the following sections, references are 
arranged in the first author's alphabetical order. Each reference 
is listed only once in the report. Effort has been made to provide 
as complete list of references as possible. 
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